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GENERAL  DISCUSSION 


In  the  following  discussion  we  will  describe  the 
work  performed  by  Science  Applications,  Inc.  (SAI)  on 
contract  number  N00014-81-C-2085  ,  SAI  project  number 
1-157-18-265,  entitled  "Fluid  Dynamics-Reactive  Flow 
Modeling,"  which  had  a  technical  performance  period  of 
1  December  1980  to  30  November  1981. 

-This  work  includes  (1)  adaptations  of  the  FAST2D 
code  for  slow  flow  applications,  (2)  use  of  the  FAST2D  code 
to  model  vorticity  and  turbulence  flow  efficiently,  (3)  an 
analytical  and  computational  study  of  the  ignition  and 
quenching  of  laminar  flames  in  a  premixed  H2 ,  0,,  N2  mix¬ 
ture,  and  (4)  incorporation  of  the  latest  reactive  flow 
computer  packages  into  the  NRL  ionosphere  code..  Each  of 
these  topics  will  be  described  in  the  following  paragraphs. 

SAI  has  recently  increased  the  flexibility  and 
efficiency  of  FAST2D,  a  two  dimensional  hydrocode  developed 
by  The  Laboratory  for  Computational  Physics  at  NRL.  The 
transport  algorithm  used  for  transporting  density,  momentum 
energy,  and  chemical  species  has  been  altered  to  compute 
along  as  many  as  ten  rows  of  computational  cells  at  a  time. 
The  original  algorithm  could  only  accommodate  one  row  at 
a  time.  The  new  algorithm  was  compared  to  the  original 
using  a  cylindrical  axisymmetric  shear  flow  problem.  Since 
no  changes  were  made  in  the  differencing  scheme  itself,  the 
new  algorithm  gave  exactly  the  same  results  as  the  original 
A  factor  of  five  increase  in  computing  speed  was  obtained. 


In  preparation  for  reactive  slow  flow  modeling 
FASTlD's  mechanism  for  following  chemical  species  was 
altered.  Number  densities  rather  than  fractional  densities 
for  all  the  species  are  transported  in  order  to  make  FAST2D 
compatible  with  existing  chemical  reaction  schemes  and 
molecular  diffusion  packages.  The  revised  version  has  been 
benchmarked  in  a  non-reacting  shear-flow  problem  similar 
to  the  one  mentioned  above.  As  a  result,  no  significant 
numerical  diffusion  of  chemical  species  occurred. 

A  project  involving  shock  modeling  was  also  completed 
using  FAST2D.  The  object  of  the  study  was  to  accurately 
and  efficiently  capture  the  various  structures  in  complex, 
non-steady  shocks.  The  quasi -Lagrangian  rezoning  capabil¬ 
ities  of  FAST2D  were  sufficient  for  modeling  two  of  the 
four  cases  studied,  regular  and  single  mach  reflection. 
Complex  and  double  mach  reflections  were  modeled  after  an 
additional  rezoning  capability  was  added  to  FAST2D.  The 
code  was  modified  to  allow  the  grid  to  expand  about  any 
arbitrary  point  within  the  computational  system.  In  this 
way,  grid  and  fluid  motions  were  synchronized  such  that 
self-similar  shock  structures  could  grow  naturally  with 
the  system.  This  work  has  been  reported  in  Fry,  M.,  et  al., 
"Shock  Capturing  Using  Flux-Corrected  Transport  Algorithms 
with  Adaptive  Gridding,"  NRL  Memo  Report  4629  (1981)  and 
appears  in  this  report  as  Appendix  A.  This  new  adaptive 
gridding  technique  allows  one  to  more  efficiently  model 
vorticity  and  turbulence  flows. 

A  theoretical  and  computational  study  of  the  igni¬ 
tion  and  quenching  of  laminar  flames  in  a  premixed  H2:02:N2 
mixture  has  been  performed  and  a  detailed  report  of  the 
work,  entitled  "A  Theoretical  Study  of  the  Ignition  of 
Pre-Mixed  Gases,"  is  included  as  Appendix  B.  The  objectives 
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of  this  work  were:  (1)  to  complete  the  calibration  of  the 
similarity  solution  model,  and  (2)  to  utilize  the  similarity 
solution  model  in  conjunction  with  the  NRL  1-D  flame  model 
to  study  the  ignition,  propagation  and  quenching  of  laminar 
flames  in  premixed  gases. 

The  similarity  solution  model  is  based  on  an 
analytic  similarity  solution  to  the  non-linear  time- 
dependent  slow  flow  equations.  The  similarity  solution 
and  the  induction  time  for  the  fuel -oxidi zer  mixture  as  a 
function  of  temperature  and  pressure  can  be  used  to  calcu¬ 
late  whether  or  not  a  given  energy  source  is  adequate  to 
ignite  the  system.  This  simple  procedure  is  then  calibrated 
using  the  N'RL  1-D  flame  model  which  includes  the  thermo¬ 
physical  properties  of  the  mixture,  a  full  chemical 
kinetics  scheme,  nonlinear  convection,  molecular  diffusion 
and  thermal  conduction.  The  details  of  the  calibration 
are  given  in  Appendix  B. 

The  similarity  solution  model  predicts  whether  or 
not  a  mixture  will  ignite  given  the  initial  radius  of  energy 
deposition  R^ ,  the  duration  of  the  heating  and  the  total 
energy  deposited  in  the  system  E  .  If  ignition  is  predicted, 
the  model  gives  the  time  it  takes  for  ignition  to  occur. 

In  contrast,  the  NRL  flame  model  not  only  predicts  ignition, 
but  also  provides  the  structure  of  the  propagating  flame. 

The  calculations  presented  in  Appendix  B  show  that  the 
similarity  solution  predicts  ignition  accurately  when  the 
radius  of  energy  deposition  is  larger  than  the  quench 
radius.  The  NRL  flame  model  is  used  to  estimate  the  quench 
radius  as  well  as  to  explain  the  discrepancy  in  the  predictions 
of  the  two  models  for  very  small  radii  of  energy  deposition. 

SAI  has  incorporated  the  latest  reactive  flow 
computer  packages  into  the  NRL  computer  code  for  the 
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numerical  modeling  of  the  mid-latitude  ionosphere . 

There  are  two  significant  additions  to  the  code  to  improve 
speed  and  flexibility.  Firstly,  the  ordinary  differential 
equation  (ODE)  solver,  IMPLCT,  used  to  time  advance  the 
chemistry  part  of  the  continuity  equations  was  replaced 
by  one  which  takes  advantage  of  the  form  of  the  equations. 

If  there  are  N'Z  grid  points  in  altitude,  then  the  ODE's 
decouple  into  N2  systems  since  with  the  time-splitting 
algorithm  used  the  chemistry  at  a  given  altitude  will  not 
depend  on  the  chemistry  at  any  other  altitude  within  the 
one  hydrodynamic  time  step.  Doing  many  systems  at  once  is 
an  ideal  situation  for  a  vector  computer  like  the  TI-ASC. 

The  computer  package,  VSAIM,  developed  at  NRL  does  such  a 
problem.  Since  the  chemistry  at  one  altitude  may  proceed 
faster  than  at  another  altitude,  VSAIM  allows  each  system 
to  use  its  own  time  step.  After  each  chemistry  time  step 
it  terminates  those  that  have  reached  the  desired  ending 
time  and  if  necessary  starts  a  new  system  thus  keeping  the 
vectors  as  optimally  long  as  possible  throughout  the  computa- 
tion.  VSAIM,  like  its  predecessor,  CHEMEQ,^"-^  also  tries 
to  treat  those  equations  it  deems  are  difficult  to  solve 
with  a  special  algorithm  for  stiff  equations  and  those 
which  are  easy  to  solve  with  a  simpler  algorithm. 

The  second  major  incorporation  into  the  ID  ionosphere 
code  is  a  routine  for  the  processing  of  chemical  reactions 
and  their  reaction  rates.  Formerly  these  reactions  and 
reaction  rates  were  hard-wired  with  the  program,  requiring 
program  modification  fpr  even  the  slightest  change  in 
reaction  rate  coefficients  or  the  number  of  reactions.  By 
using  the  automatic  rate  processor,  ARTP,  the  chemical 
reactions  and  their  corresponding  reaction  rates  are  input 
as  data  cards  allowing  great  flexibility  in  the  choice  of 
the  type  of  chemistry  one  elects  to  incorporate  into  the 
ionospheric  model. 
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In  order  to  make  these  modifications  it  was  neces¬ 
sary  to  change  the  data  structures  used  to  store  the  densities 
of  the  various  chemical  species.  Instead  of  using  a 
distinct  array  for  each  specie  present  in  the  problem, 
one  larger  array  indexed  by  the  number  of  species  in  the  ■ 
problem  is  now  used.  An  array  of  symbols  is  now  used  to 
distinguish  one  specie  from  another.  These  symbols  are 
the  data  that  ARTP  uses  to  recognize  which  species  are 
involved  in  each  reaction  it  processes.  By  making  the 
types  of  ions  and  minor  neutrals  input  to  the  program, 
the  user  has  a  greater  flexibility  to  test  the  effects  of 
various  chemical  reactions  and  the  importance  of  the  various 
ions  or  minor  neutrals  on  the  resulting  ionosphere. 

While  the  modifications  to  the  code  were  being  made, 
documentation  describing  the  variables  used  by  the  code  was 
added.  Furthermore,  a  separate  piece  of  documentation  was 
added  to  describe  the  input  data  to  the  program,  making 
the  whole  program  much  more  user  oriented.  In  .Appendix  C 
we  have  included  the  revised  versions  of  2  subroutines  -  I.XITAL 
and  RTCO.N.  These  subroutines  supply  all  the  information 
necessary  for  setting  up  the  input  required  for  making  a 
run  of  the  ionosphere  code. 
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A  numerical  technique  has  been  developed  for  capturing  com¬ 
plex,  nonsteady  shock  structures  in  multidimensions.  The 
technique  relies  on  moving  the  computational  mesh  with  the 
shock  wave  so  chat  Che  features  of  principal  interest  appear 
approximately  stationary.  The  method  has  been  implemented 
using  coordinate-split  Flux-Corrected  Transport  (FCT)  al¬ 
gorithms  which  allow  the  mesh  to  evolve  arbitrarily  with 
respect  Co  the  fluid  in  each  coordinate.  The  grid  may  thus 
be  optimized  in  response  to  the  needs  of  a  given  problem. 
Synchronizing  Che  grid  and  fluid  motions  permits  significant 
reduction  of  numerical  transients  and  eliminates  numerical 
diffusion.  Shocks  develop  naturally,  with  no  fitting.  The 
method  is  illustrated  by  calculating  complex,  two-dimension¬ 
al  Mach  reflection  phenomena  associated  with  airblasts  and 
shock  diffraction  on  wedges.  The  numerical  results  are  in 
good  agreement  with  available  experimental  data. 

INTRODUCTION 

Numerical  solution  of  transient  multidimensional  gas  dynamics  problems  is 
always  nontrivial.  When,  in  addition,  the  problem  involves  reflecting  super¬ 
sonic  flows,  large  variations  in  length  scales  in  both  space  and  time,  or  phe¬ 
nomena  for  which  neither  analytic  solutions  nor  detailed  experimental  observa¬ 
tions  are  at  hand,  the  state  of  the  computational  art  is  challenged.  Such  a 
problem  arises  in  calculating  the  oblique  reflection  of  shocks  from  solid  sur¬ 
faces  in  planar  geometries  (e.g.  shock  tube  experiments)  or  axisymmetric  geo¬ 
metries  (e.g.  airblasts).  The  complications  arise  mainly  from  the  presence 
of  Mach  reflections  which  occur  when  a  shock  front  impinges  on  a  reflecting 
surface  at  angles  of  incidence  sufficiently  far  from  normal.  The  formation  of 
a  Mach  stem  and,  consequently,  of  a  slip  surface  intersecting  the  triple  point 
(the  confluence  of  the  Incident,  Mach,  and  reflected  waves)  results  from  the 
requirement  that  the  flow  behind  the  reflected  shock  be  parallel  to  the  re¬ 
flecting  surface,  which  cannot  be  achieved  through  regular  reflection. 


explosive  studies^'^.  For  the  simplest  problem,  that  of  a  planar  shock 


Manuscript  aubmitted  August  24, 1981. 


A- 2 


t 


( 


reflecting  from  a  plane  surface,  Jones,  Martin,  and  Thornhill^  noted  that  it 
is  possible  to  reduce  the  number  of  independent  variables  to  two  by  transform¬ 
ing  to  the  similarity  variables  xjt,  y/t,  a  device  that  was  also  used  by 
Kutler,  et  al^.  Ben-Dor^  developed  a  theory  which  used  shock  polars  to  explain 
some  of  the  features  of  this  problem,  and  solved  the  system  of  algebraic  equa¬ 
tions  obtained  by  combining  the  Jump  conditions  across  the  various  disconti¬ 
nuities  (Courant  and  Friedrichs)”  to  describe  the  flow  in  the  neighborhood  of 
the  triple  point.  To  date,  no  satisfactory  treatment  of  the  complete  flow 
field  has  been  published,  although  some  features  (like  the  shape  of  various 
waveforms)  are  quite  easy  to  model. 


In  connection  with  studies  of  both  chemical  and  nuclear  explosions  there 
have  been  many  attempts  to  model  a  spherical  blast  wave  reflecting  from  the 
ground,  the  so-called  height-of-burst  (HOB)  problem.  The  hydrodynamic  pheno¬ 
mena  in  the  two  cases  are  identical,  although  nonideal  effects  (primarily  ex¬ 
plosive  afterburn  in  the  first  instance  and  radiation  preheating  in  the  second) 
are  different.  Previous  attempts  to  model  two-dimensional  complex  shock  re¬ 
flection  have  suffered  from  restriction  to  describing  part  of  the  system,  the 
use  of  a  special  assumption  like  that  of  self-similarity,  or  less  than  satis¬ 
factory  agreement  with  experimental  data.^ 


The  calculations  discussed  here  represent  a  step  forward  in  overcoming 
these  difficulties.  They  differ  from  previous  nximerical  work  in  incorporating 
two  important  computational  developments:  Flux-Corrected  Transport  (FCT)^^  and 
an  adaptive  regrlddlng  procedure,  called  "sliding  rezone", which  optimizes 
the  mesh  point  distribution  and  hence  the  resolution  of  surfaces  of  disconti¬ 
nuity. 


FCT  is  a  finite-difference  technique  for  solving  the  fluid  equations  in 
problems  where  sharp  discontinuities  arise  (e.g.  shocks,  slip  surfaces  and 
contact  surfaces) .  It  modifies  the  linear  properties  of  a  second-  (or  higher) 
order  algorithm  by  adding  a  diffusion  term  during  convective  transport,  and 
then  subtracting  it  out  "almost  everywhere"  in  the  antidiffusion  phase  of  each 
time  step.  The  residual  diffusion  is  just  large  enough  to  prevent  dispersive 
ripples  from  arising  at  the  discontinuity,  thus  ensuring  that  all  conserved 
quantities  remain  positive.  FCT  captures  shocks  accurately  over  a  wide  range 
of  parameters.  Ho  information  about  the  number  or  nature  of  the  surfaces  of 
discontinuity  need  be  provided  prior  to  initiating  the  calculation. 


The  FCT  routine  used  in  the  present  calculations,  called  JPBFCT  (an  ad¬ 
vanced  version  of  ETBFCT) consists  of  a  flexible,  general  transport  module 
which  solves  1-D  fluid  equations  in  Cartesian,  cylindrical,  or  spherical  geo¬ 
metry.  It  provides  a  finite  difference  approximation  to  the  conservation  laws 
of  the  general  form:. 


0  f  MV  “  -/*  ♦  (u-u  )  ■  dA  +  y  TdA  (1) 

dt  76V(t)  J  bA(t)  ''«A(t) 

where  represents  the  mass,  momentum,  energy  or  mass  species  in  cell  6V(t), 

II  and  ^  represent  the  fluid  and  grid  velocities,  respectively,  and  t  repre¬ 
sents  ^e  pressure/work  terms.  This  formulation  allows  the  grid  to  slide  with 
respect  to  the  fluid  without  introducing  any  additional  numerical  diffusion. 
Thus,  knowing  where  the  features  of  greatest  Interest  are  located,  one  can 
concentrate  fine  zones  where  they  will  resolve  these  features  most  effectively 
as  the  system  evolves  (Fig.  1). 

In  the  next  section  we  describe  the  computational  techniques  used  to  solve 
the  wedge  problem  and  present  the  results  of  four  simulations  carried  out  to 
reproduce  experimental  results  of  Ben-Dor  and  Glass. In  Section  III  we  pre¬ 
sent  a  parallel  discussion  for  a  HOB  calculation.  Finally,  in  Section  IV  we 
summarize  our  conclusions. 


Fig.  1.  Adaptive  grids  for  a)  planar  shocks  on  wedge  (double  Mach  shock  fea¬ 
tures  are  Indicated) ;  b)  and  c)  HOB  problem  Initially  and  at  transition  point 
(grid  lines  In  fine-zone  region  are  Indistinguishable) . 

SHOCK-ON-WEDGE  CALCULATIONS 

The  JPBFCT  algorithm  was  used  In  a  2-D  Cartesian  version  of  the  FAST2D  code 
to  model  the  reflections  of  planar  shocks  from  wedges  of  20<>  to  60°  and  vary¬ 
ing  shock  strengths.  Four  general  classes  which  include  regular,  single,  com¬ 
plex  and  double  Mach  reflection  were  calculated  (referred  to  as  cases  a,b,c,d 
respectively) .  The  bottom  of  the  mesh,  treated  as  a  reflecting  boundary, 
modeled  the  surface  of  the  wedge.  Quantities  on  the  right  hand  boimdary  and 
on  the  top  were  set  equal  to  the  ambient  values.  The  remaining  boundaries  were 
treated  as  permeable.  In  the  single,  complex,  and  double  Mach  reflection 
cases,  the  mesh  was  anchored  on  the  left,  essentially  at  the  wedge  tip  where 
the  incident  shock  first  strikes,  while  the  zones  were  stretched  by  a  scaling 
factor  proportional  to  t  as  soon  as  the  reflection  region  filled  a  substantial 
portion  of  the  grid.  In  case  (d),  the  double  Mach  reflection  case,  the  open¬ 
ing  angle  is  so  small  that  the  Incident  shock  has  to  traverse  many  zones 
before  the  mach  stem  has  grown  large  enough  to  be  well  resolved.  For  this  rea¬ 
son,  the  problem  was  solved  on  a  uniform  mesh  in  the  frame  of  reference  fixed 


to  the  reflection  point,  with  stretching  being  initiated  after  the  first  Mach 
stem  reached  20  cells  in  length.  The  timestep  was  recalculated  at  every 
cycle  with  a  Courant  number  of  0.5. 

Figure  2  shows  the  pressure  and  density  contours  and  the  velocity  field  for 
cases  a,b,c,d.  The  pertinent  shock  phenomena  can  be  easily  identified:  inci¬ 
dent  shock,  contact  surface,  first  and  second  Mach  stems.  As  shown  in  Fig.  1, 
the  zoning  is  particularly  sparse  except  for  the  region  of  Interest.  Adequate 
resolution  of  the  key  surfaces  (contact  and  second  Mach  stem)  is  obtained  with 
5  zones  in  each  direction.  The  accuracy  can  be  evaluated  by  comparing  the  ex¬ 
perimental  density  distributions  along  the  wall  (Fig.  3). 
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Fig.  3.  Comparison  of  density  (In  units  of  ambient  density  p  )  for  cases  (a), 
(b),  (c),  (d)  of  Fig.  2  vs.  distance  from  comer.  Points  are  measured  values 
reported  In  Ref.  13. 

HEIGHT  OF  BURST  CALCULATIONS 

Next,  we  performed  a  numerical  simulation  of  a  IKT  nuclear  detonation  at 
31.7  m  HOB,  a  case  which  could  be  readily  compared  with  high  explosive  data.  A 
constant  ambient  atmosphere  was  used  with  a  density  of  1.22  x  10~3  g/cm^  and 
pressure  1.01  x  10®  dynes/cm^.  To  relate  the  energy  and  density  to  the  pressure, 
a  real-air  equation  of  state  (EOS)  was  used.  This  table-lookup  EOS  was  derived 
from  theoretical  calculations  by  Gllmore^^’^^  for  equilibrium  properties  of  air 
and  has  been  vectorized  for  the  Advanced  Scientific  Computer^®.  The  internal 
energy  density  used  In  the  call  to  the  EOS  Is  found  by  subtracting  kinetic 
from  total  energy;  this  can  be  negative  due  to  truncation  (phase)  errors.  When 
this  occurred,  the  value  of  the  pressure  was  reset  to  zero. 

The  transition  from  regular  reflection  to  double  Mach  reflection  occurs 
at  a  ground  range  approximately  equal  to  the  HOB.  The  size  of  the  mesh  should 
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therefore  be  roughly  twice  the  HOB  in  both  directions.  The  upper  boundary 
should  be  far  enough  away  from  the  blast  front  to  be  non-lnterferlng.  We  chose 
boundaries  of  55  m  for  the  radial  direction  and  103.5  m  for  the  axial  direction. 
The  fine  grid  In  the  radial  direction  contained  140  out  of  200  total  zones, 
each  5  cm  In  length.  The  rightmost  zones  were  80  cm  In  length,  and  a  smoothing 
Involving  40  zones  performed  between  the  regions  to  guarantee  that  the  zone 
sizes  varied  slowly.  In  the  axial  direction  the  fine  grid  contained  75  out  of 
150  total  zones,  each  5  cm  In  length.  Beyond  that  region  the  zones  were  geo¬ 
metrically  increased  by  a  factor  of  1.112. 

Placement  of  the  fine  grid  at  the  origin  of  the  mesh  (ground  zero,  the 
point  at  which  reflection  first  occurs)  was  determined  to  be  optimum  for  cap¬ 
turing  peak  pressure  In  the  alrblast  wavefront.  Thus,  as  the  expanding  wave 
moves  along  the  ground  surface,  the  fine  grid  Is  always  locked  to  It  and  each 
point  along  the  blast  front  encounters  the  same  spatial  grlddlng  as  It  approaches 
the  ground.  By  treating  each  point  of  the  Incident  front  In  the  same  manner, 
we  Insure  that  the  calculation  Is  Internally  consistent  and  that  the  computed 
transition  point  Is  accurate  to  within  the  limits  of  the  resolution. 

The  Initialization  provides  a  strong  shock  with  approximate  Mach  number 
M»12.  This  speed  and  the  need  for  restart  capability  led  to  the  choice  of  200 
tlmesteps  as  an  interval  for  the  spatial  display  (snapshots) .  The  dump  Inter¬ 
val  that  resulted  was  At  0.3  milliseconds  (ms).  These  dumps  were  stored  on 
magnetic  tape  and  post-processed. 

A  fit  to  the  1-D  nuclear  blast  flow  field  (Ref.  17)  was  used  to  Initialize 
the  energy  and  mass  density  and  velocity  field  at  3.76  ms.  The  corresponding 
peak  overpressure  was  113  bars.  After  the  1  KT  flow  field  was  laid  down  Inside 
a  radius  of  31.6  m,  the  fine-zone  grid  was  activated  to  follow  the  peak  pressure 
as  It  moved  along  the  ground  surface,  modelled  as  a  perfectly  reflecting  bound¬ 
ary.  This  region  comprised  140  zones,  and  a  switch  was  set  to  keep  40  of  these 
zones  ahead  of  the  reflection  point.  Permeable  boundary  conditions  are  used 
on  the  top  and  right  edges  of  the  mesh,  i.e.,  density,  pressure  and  velocity 
are  set  equal  to  ambient  preshock  conditions.  Reflecting  conditions  were 
applied  to  the  left  and  bottom.  The  total  elapsed  physical  time  in  the  2-D 
calculation,  7.6  ms,  required  5600  cycles.  Times  are  referred  to  t*0  at  the 
start  of  the  calculation. 


The  numerical  simulation  begins  just  before  the  shock  first  reflects  from 
the  ground.  Fig.  4a  Indicates  the  pressure  and  density  contours  and  velocity 
vectors  at  time  3.18  ms.  In  Fig.  4b  the  reflected  shock  Is  shown  moving  upward, 
the  outward  flow  begins  to  stagnate  at  the  ground  (transition).  Fig.  4c,  t*5.99 
ms,  shows  an  enlargement  of  the  shockfront,  and  the  development  of  the  Mach 
stem,  slip  surface  and  second  Mach  stem.  The  angle  of  the  shock  front  with 
respect  to  the  ground  Is  Increasing  with  time,  so  that  the  effective  wedge 
angle  is  decreasing.  From  Ben-Dor  and  Glass^^  one  expects  a  transition  to 
double  Mach  stem  to  occur  at  approximately  45**.  The  angle  In  Fig.  4b  Is  about 
45*^  and  the  shock  front  has  entered  the  transition  phase.  Figure  4d  shows  the 
fully  developed  shock  structure  at  7.79  ms.  Clearly  visible  is  the  second  Mach 
stem  and  a  vortex  region  behind  the  first  Mach  stem.  Toeing  out  of  the  first 
Mach  stem  can  be  also  seen  In  the  contours  of  Fig.  4d  and  occurs  as  the  fluid 
rolls  forward  where  the  slip  line  would  otherwise  Intersect  the  ground.  The 
velocity  field  In  Fig.  4d  also  shows  this  detail. 

One  should  also  note  the  reflected  shock  properties.  The  reflected  shock 
propagates  rapidly  through  the  high  temperature  fireball,  due  to  the  high  local 
soundspeed.  The  shape  of  this  reflected  wave  Is  a  primary  difference  between 
the  HOB  case  and  the  wedge  case^^.  The  other  major  difference,  of  course.  Is 
the  spherically  expanding  blast  wave  which  decreases  In  strength  approximately 
proportional  to  r~^. 


PRESSURE 


density 


VELOCITY 


Fig.  4.  Pressure,  density,  and  velocity  fields  for  HOB  calculation  (a)  in 
regular  reflection  stage;  (b)  at  transition  to  Mach  reflection;  (c)  shortly 
afterward,  when  second  peak  has  become  larger  than  first;  and  (d)  fully  de¬ 
veloped  (note  toe  at  base  of  first  Mach  stem) . 


Finally  we  consider  the  pressure/distance  relation  for  the  HOB  case.  In 
Fig.  5  we  compare  the  results  of  the  numerical  simulation  with  the  data  of 
Carpenter  and  with  empirical  analysis.  Carpenter's  data  are  based  upon  care¬ 
ful  HOB  experiments  with  8  lb  FBX9404  spheres.  The  empirical  analysis  was 
based  on  a  1  KT  nuclear  free  air  curve  and  HOB  construction  factors.  The  cal¬ 
culated  values  in  the  regular  reflection  regime  are  20Z  low  and  may  be  attri¬ 
buted  to  a  combination  of  FCT  clipping,  the  resolution  of  the  grid,  and  in¬ 
accuracies  in  the  initialization  of  the  flow  field.  During  and  after  Mach 
reflection,  the  peaks  remain  low  until  the  Mach  stem  structure  has  grown  large 
enough  to  be  resolved  on  the  mesh.  By  the  time  it  occupies  a  region  of  15  cells 
high  and  35  cells  wide,  the  peak  pressures  are  in  good  agreement  with  the  HE 
data  and  the  empirical  analysis. 
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Fig.  5.  Pressure-range  curves  for  first  and  second  (after  transition  -  denoted 
by  TP  -  to  double  Mach  reflection)  peaks. 

SUMMARY  AND  CONCLUSION 

The  complex  2-D  Mach  reflection  phenomena  associated  with  shock  diffraction 
on  wedges  and  helght-of-burst  explosions  have  been  modeled  with  the  FAST2D  com¬ 
puter  code.  Four  wedge  cases — regular,  single,  complex  and  double  Mach  reflec¬ 
tion — have  been  calculated  aixi  the  results  compared  to  experiments.  A  nuclear 
detonation  (1  KT  at  31.7m  HOB)  was  also  simulated.  The  results  give  insight  in¬ 
to  the  formation  and  subsequent  evolution  of  the  Mach  stem,  the  triple  point 
and  the  contact  discontinuity.  The  transition  from  regular  reflection  to  double 
Mach  reflection  is  predicted.  Excellent  agreement  with  Ben-Dor's  data  is  obtained. 
We  suggest  that  the  first  signal  for  transition  Is  the  appearance  of  a  second 
peak  behind  the  shock  front  due  to  stagnation  in  the  flow.  Calculated  first  and 
second  pressure  peaks  versus  distance  in  the  HOB  case  agree  both  with  the  HE 
data  and  analysis  to  within  20Z. 

The  use  of  the  adaptive  regridding  procedure,  called  "sliding  rezone", along 
with  the  FCT  algorithm  allows  one  to  accurately  predict  the  nonsteady  shock 
structures  in  two  dimensions  for  diffractions  on  wedges  and  HOB  cases.  Compari¬ 
son  with  data  for  both  wedges  and  HOB  yields  the  best  results  obtained  to  date. 
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ABSTRACT 


In  this  paper,  time-dependent  results  obtained  from  both  a  simple 
but  nonlinear  analytic  similarity  solution  and  a  detailed  numerical  simula¬ 
tion  model  are  used  to  study  the  interactions  between  the  fundamental  processes 
occurring  in  the  ignition  of  homogeneous  premixed  geises.  The  parameters  which 
may  be  varied  are  the  composition  of  the  mixture,  the  initial  radius  of  energy 
deposition  R^,  the  duration  of  the  heating  and  the  total  energy  deposited 
in  the  system  E^.  The  similarity  solution  plus  the  Induction  time  for  the 
fuel-oxidizer  mixture  as  a  function  of  temperature  and  pressure  can  be  used 
to  calculate  whether  or  not  a  given  energy  source  is  adequate  to  ignite  the 
system.  This  simple  procedure  is  then  calibrated  using  a  time-dependent 
detailed  numerical  reactive  flow  model  which  includes  the  thermophysical 
properties  of  the  mixture,  a  full  chemical  kinetics  scheme,  nonlinear  convec¬ 
tion,  molecular  diffusion  and  thermal  conduction.  Calculations  are  presented 

for  a  selected  mixture  of  H_-0.,-N.,  for  various  values  of  R  and  E  .  These 

2  2  2  O  O 

show  that  the  similarity  solution  predicts  ignition  accurately  when  the  radius 
of  energy  deposition  is  larger  than  the  quench-radius .  The  detailed  numerical 
reactive  flow  model  is  used  to  predict  the  quench-radius  and  the  absolute 
minimum  ignition  energy  associated  with  it. 
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Nonlinear  amplitude  in  similaurity  solution 
Energy 

Total  amount  of  energy  deposited 

Induction  Parameter 

Scale  size  in  similarity  solution 

Boltzmann's  const2uit 

Total  number  density 

Number  density  of  species  j 

Pressure 

Term  representing  production  of  species  j 

Contribution  to  the  heat  flux  vector  due  to  diffusion 

Term  representing  loss  of  species  j 

Initial  radius  of  deposition 

Characteristic  radius  of  energy  deposition 

Radial  distemce 

Tenqperature 

time 

fluid  velocity 

Diffusion  velocity  of  species  j 
ratio  of  specific  heats,  c  /c 

p  V 

mixture  viscosity  coefficient 

mixture  thermal  conductivity  coefficient 

mass  density 

time  period  during  which  energy  is  deposited 


chemical- induction  time 


introduction 


An  external  source  of  energy  cam  initiate  interactions  among  the 
controlling  convective,  diffusive  and  chemical  processes  In  a  fuel- 
oxidizer  mixture.  Whether  the  interactions  result  in  ignition  of  a 
deflagration  or  detonation  wave  depends  on  the  intensity,  duration, 
and  volvime  affected  by  the  external  heat  source.  Ignition  will  also 
depend  on  the  initial  ambient  properties  of  the  mixture  which  determine 
the  chemical  induction  time  and  the  heat  release  per  gram  of  material. 

Thus  ignition  is  a  complicated  phenomena  whose  occurrence  for  a  specific 
mixture  of  fuel  and  oxidizer  depends  strongly  on  diffusive  and  chemical 
parcuneters  which  are  often  very  poorly  known. 

It  is  possible,  in  principle,  to  study  flame  ignition  by  performing 
detailed  numerical  simulations.  This  is  a  complicated,  multi-dimensional, 
multi-species,  time-dependent  problem  which  has  been  solved  in  certain  limited 
cases  (see  e.g..  Ref.  1).  Part  of  the  complication  and  cost  of  such  cal¬ 
culations  arises  from  the  multi-dimensional  solution  of  the  conservation 
equations ,  but  at  least  as  much  arises  from  integrating  the  large  number  of 
ordinary  differential  equations  describing  the  chemical  reactions.  This 
latter  factor  is  further  complicated  by  the  fact  that  we  usually  do  not  have 
an  adequate  representation  of  the  chemical  reactions  with  which  to  work.  Thus 
a  convenient,  inexpensive  way  to  estimate  whether  a  mixture  will  ignite  given 
a  heat  source  intensity,  duration,  and  volume  would  be  a  very  valuable  laborer 
tory  tool  and  a  useful  learning  device. 
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This  paper  represents  an  extension  and  clarification  of  the  work 
described  by  Oran  and  Boris  [1,2].  They  presented  a  preliminary  summary 
of  a  simplified,  theoretical  model  of  localized  ignition  of  a  homogeneous, 
premixed  gas  and  explained  how  it  could  be  calibrated  by  using  a  detailed 
simulation.  In  this  paper  the  calibration  of  the  simplified  model  is 
extended  and  the  results  are  conpared  in  detail  to  those  of  detailed, 
one-dimensional  simulations.  After  a  brief  review  of  related  literature, 
the  essential  properties  of  the  detailed,  time-dependent,  numerical  model 
are  presented.  Then  a  description  is  given  of  a  closed  form  similarity 
solution  for  the  nonlinear  time-dependent  slow  flow  equation  which  forms  the 
basis  for  the  simplified  model  of  localized  ignition.  The  similarity  model 
avoids  the  integration  of  the  ordinary  differential  equations  describing 
the  chemical  reactions  by  defining  an  induction  parameter.  Two  constants 
must  be  calibrated;  the  radius  at  which  the  thermal  conductivity  is  evaluated 
and  the  radius  at  which  the  induction  parameter  is  evaluated.  Finally, 
comparisons  between  the  two  models  are  described  for  ignition  of  a  H2“02“^2 
mixture.  The  composition  of  the  mixture  and  the  duration  of  energy  deposition 
are  held  constant  but  the  total  energy  deposited  and  the  radius  of  deposition 


are  allowed  to  vary. 


BACKGROUND 


Ignition  phenomena  ^md  the  associated  properties  of  minimum  ignition 
energies  and  quench  volumes  have  been  studied  both  experimentally  and 
theoretically.  Lewis  and  von  Elbe  [3]  have  reported  extensive  experimental 
data  on  ignition  by  electric  sparks.  Weinberg  etnd  Wilson  [4]  and  Kingdon 
and  Weinberg  [5]  have  compared  the  ignition  energies  required  when  a  spark 
and  a  Iciser  were  the  source  of  ignition.  The  former  paper  [4]  concludes  that 
the  laser  minimum  ignition  energy  is  very  much  less  than  that  required  by  a 
spark  for  mixtures  at  low  pressures  or  near  the  flammability  limit.  They  have 
attributed  this  difference  to  the  influence  of  the  electrodes  and  the  losses 
to  them  since  the  quenching  distemces  are  large  under  the  adsove  conditions. 

The  latter  paper  [5]  concludes  that  for  short  pulses  applied  to  mixtures  with 
small  quenching  distances,  the  energies  required  by  the  two  ignition  sources 
are  not  that  different. 

Dixon-Lewis  and  Shepherd  [6]  have  used  a  time-dependent  flame  model  to 
examine  the  effects  of  varying  the  initial  radical  concentration  and  of 
changing  the  geometry  of  the  initial  hot  core.  Their  studies  were  done  on  a 
homogeneous  premixed  60%  hydrogen-air  mixture  in  which  energy  was  deposited 
Instantaneously  in  the  form  of  hydrogen  atom  radicals  and  temperature.  They 
found  that  energy  in  the  form  of  hydrogen  atoms  was  more  efficient  at  igniting 
a  flcune.  Dixon-Lewis  [7]  then  looked  for  a  minimum  quench  volume  in  the  same 
60%  hydrogen-air  mixture  in  which  one  third  of  the  added  energy  was  in  the 
form  of  hydrogen  atoms. 
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Another  point  noted  by  both  Over ley  et  el.  [8]  and  by  Dixon-Lewls 


and  Shepherd  [6]  is  that  the  shape  and  size  of  the  mixture  core  is 
important  in  deteannining  the  minimum  ignition  energy.  The  former  paper 
presented  a  study  of  a  hydrazine  mixture  in  both  spherical  euid  prolate  . 
spheroidal  coordinates.  The  latter  paper  considered  both  cylindrical  £uid 
spherical  systems.  Both  papers  investigated  the  effects  of  flame  curvature 
on  burning  velocity. 

The  work  of  Ballal  and  Lefebvre  [9]  is  a  detailed  eacperimental  investi¬ 
gation  of  the  effects  of  pressure,  velocity,  turbulent  intensity  and  scale, 
and  mixture  composition  on  the  minimum  ignition  energy  and  quench  distance 
in  a  flowing  gas  mixture.  Using  the  data  as  a  guide,  they  analyzed  the 
important  transport  and  turbulent  processes  and  developed  a  highly  simplified 
model  for  the  quench  distance  in  the  lim..c  of  low  and  high  turbulence.  Hie 
data  has  allowed  them  to  calibrate  the  model  and  determine  important  constants. 
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THE  DETATT-BD  NUMERICAL  SIMULATION  MODEL 

This  section  describes  the  detailed  numerical  flame  model  which  solves 
the  time-dependent  conservation  equations  for  total  mass  p,  momentum  pv» 
energy  E,  euid  individual  species  densities  n^  [2,10]: 

1^  .  -V  .  pv  (1) 

3n. 

-r-i  =  -7  •  n.V.  -  7  •  n.v  +  P.  -  Q.n.  (2) 

at  -  -  D-  3  *^1  3 

3Pv  _ 

— —  =  -V  •  (pw)  -  7P  +  7  •  n„[7v  +  (7v)  1  (3) 

II  =  -7  •  Ev  -  7  .  (Pv  -  Ijj^T  -  ^)  .  (4) 

The  technique  for  solving  the  various  terms  in  these  equations  is  based 
on  the  method  of  asymptotic  timestep  splitting  in  which  the  individual 
chemical  and  physical  processes  are  integrated  sepaurately  and  then  asymptot¬ 
ically  coupled  together  [1],  The  model  permits  a  wide  variety  of  geometric, 
initial,  boundary  and  time-varying  energy  input  conditions  and  was  specif¬ 
ically  developed  to  study  the  various  physical  and  chemical  processes  which 
control  flame  initiation  and  quenching. [23]. 

T^e  convective  transport  is  solved  by  the  algorithm  ADINC,  a  Lagrangian 
hydrodynamic  algorithm  which  solves  implicitly  for  the  pressures  [11].  Ihe 
method  allows  for  arbitrary  equations  of  state,  gives  an  accurate  representa¬ 
tion  of  material  interfaces,  amd  allows  steep  gradients  in  species  and 
temperature  to  be  developed  and  maintained. 
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A  nvunber  of  tests  of  this  algorithm  have  been  dociunented  by  Boris  [11]  •  An 
adaptive  gridding  method  h£is  been  implemented  in  which  cells  are  inserted  or 
deleted  according  to  externsLlly  specified  physical  conditions  in  the  flow. 

The  diffusion  velocities  are  calculated  using  an  iterative  algorithm  [12, IJ. 

The  chemical  interactions  aure  described  by  a  set  of  nonlinear  coupled 
ordinaury  differential  equations.  For  this  ignition  study  we  have  used  the 

reaction  scheme  [13]  given  in  Table  1  which  involves  the  eight  reactive 
species,  O^,  O,  H,  OH,  HO^,  ^2^2’  diluent  which  is  chosen  to  be 

Nj.  The  thermochemical  properties  of  these  species  were  taken  from  the  JAMAF 
tables  [14].  The  ordinary  differential  equations  describing  the  chemical 
kinetics  are  solved  tising  a  fully  vectorized  version  of  the  selected  asyn^totic 
integration  method,  CHEMEQ,  developed  by  Young  [15,16]. 

Results  from  a  typical  flame  calculation  in  spherical  geometry  are 

presented  in  Fig.  1.  The  figure  depicts  the  time  history  of  the  temperature 

-4 

profile  after  4mJ  of  energy  is  deposited  over  a  period  of  10  seconds.  Even 
after  the  energy  deposition  is  stopped  the  "core"  temperature  continues  to 
increase  due  to  the  heat  released  in  chemical  reactions.  With  time,  however, 
the  temperature  near  the  center  decreases  and  the  temperature  away  from  the 
center  increases  due  to  heat  conduction.  In  this  case,  since  the  rate  of 
heat  generation  is  greater  than  the  rate  of  heat  loss,  the  temperature  profile 
develops  into  that  of  a  typical  flame  temperature  profile.  Species  profiles 
and  details  of  the  fl2une  front  have  been  presented  elsewhere  [2]. 


THE  SIMILARITY  SOIOTION 

The  similarity  solution  is  derived  from  the  slow  flow  approximation 
[  1,171  predicated  on  the  assumption  that  energy  addition  to  the 

system  is  slow  enough  so  that  there  is  no  shock  heating.  Thus  the  system 
is  characterized  flow  velocities  which  are  small  coii{>ared  to  the  speed 
of  sound ,  and  an  essentially  constant  pressure  field.  The  energy  and 
velocity  equations  may  be  written  in  terms  of  total  derivatives  as 

2  2 

^  AS  0  ^  -  yPV  *  V  +  7  *  +  s(t)e~^  (5) 


and  the  continuity  equation  is  written  in  the  form 


P  dt 


V. 


(6) 


In  Eq.  (5)  Y  is  held  consteuit,  and  ic  is  a  function  of  the  mixture  thermal 

conductivity,  X  , 
n 

B 

The  last  term  on  the  right  hand  side  of  Eq.  (5)  is  the  source  term.  Proper 
choice  of  S(t)  ensures  that  a  given  amount  of  energy,  E^,  is  deposited  in 
a  certain  time,  T^.  It  is  the  choice  of  this  Gaussian  profile  which  allows 
us  to  obtain  the  "closed"  form  similarity  solution  which  is  given  below. 

The  details  of  the  solution  of  Eqs.  (5)  and  (6)  are  given  in  the 


Appendix  at  the  end  of  this  paper.  The  final  results  for  the  ten^rature 
and  density  as  a  function  of  time  and  position  may  be  written  as 


t 


T(r,t) 


T^e 


A(t)e 


2  2 
-k‘"(t)r 


(8) 


and 


p (r,t) 


Pao« 


-A(t)e 


-k^(t)r^ 


(9) 


where  and  are  the  background  temperature  and  density,  respectively. 
In  addition,  two  ordinary  differential  equations  describing  the  evolution 
of  A(t)  and  k(t)  are  required  to  complete  the  solution. 


dt 


-  kv. 


2ick^ 


(10) 


M  =  S(t) 

dt  “  yP 

m 


2ck^  A 


(11) 


where  v^  is  an  approximation  to  the  velocity  based  on  the  assumption  that 
v(r,t)  =  Vj^(t)r;  and  c  is  a  constant  which  depends  upon  the  configuration 
of  the  problem  (c  =  1  for  cartesian,  2  for  cylindrical  emd  3  for  spherical 
coordinates) .  By  comparing  the  results  from  this  similarity  solution  model 
to  those  from  the  detailed  numerical  simulation  model  described  earlier,  we 
have  seen  that  the  linear  dependence  on  r  is  a  valid  approximation  before 
ignition  occurs. 

The  model  requries  one  further  defintion  in  order  to  predict  ignition. 
The  chemical  induction  time,  which  is  a  function  of  temperature  and  pressure, 
must  be  used  to  define  the  induction  parameter. 


r  •St’ 

(T,P)  -I - 

T  [T(r,f),  P(r,f 
0  ® 


(12) 
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For  the  work  presented  in  this  paper,  pressure  is  constant  throughout  the 
calculation.  Then  is  integrated  in  time  according  to  Eq.  (12). 

Ignition  "occurs"  when  I{T,P^)  ■  1,  which  is  an  exact  result  in  the 
limit  of  large  heat  source  and  constant  temperature  near  the  center  of  the 
heated  region.  Values  of  t  have  been  calculated  [13]  for  a  wide  range  of 
temperatures  imd  pressures  from  the  chemical  rate  scheme  given  in  Table  I. 

In  general,  as  required  in  these  calculations  may  be  obtained  from  such 
detailed  kinetic  calculations,  a  few  measured  points,  or  from  educated  guesses. 
The  concept  of  an  Induction  parameter  has  been  extended  in  scope  and  used 
extensively  for  reactive  shock  and  detonation  studies  [e.g.,  1,  18,  22]. 

The  parameter  k  in  Eqs.  (10)  and  (11)  is  found  using  the  definition 
in  Eq.  (7) .  Determining  the  temperature  (radius)  at  which  k  is  to  be  evaltiated 
is  part  of  the  calibration  to  be  done.  For  this,  the  derivation  in  the 
appendix  was  also  done  for  the  case  where  k  was  assumed  to  be  a  function  of 
radius  (r) .  This  analysis  indicated  that  if  k  is  to  maintain  its  similarity 
form  and  be  a  function  of  time  alone,  k.  must  be  independent  of  r.  Thus,  within 
the  framework  of  the  similarity  solution,  the  appropriate  location  at  which 
to  evaliiate  <  is  at  the  center  (r  *  0) .  Determining  the  location  at  which 
the  induction  parameter,  l(T,P^),  should  be  evaluated  is  the  second 
calibration  to  be  done.  This  has  been  done  by  comparing  the  predictions 
of  the  two  models  and  is  discussed  in  the  next  section. 

In  summary,  the  approximations  to  Eqs.  (1)  -  (4)  which  allows  us  to 
write  Eq.  (S)  and  obtain  the  solutions  represented  by  Eqs.  (8)  -  (11)  are: 
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1.  The  flow  velocities  characteristic  of  the  system  are  small  con^ared  to 
the  speed  of  sound: 

2.  The  pressure  is  essentially  constant  everywhere; 

3.  The  ratio  of  specific  heats,  y,  is  constant; 

4.  Molecular  diffusion  effects  are  not  important  until  after  ignition; 

5.  The  velocity,  V(r,t)  may  be  approximated  by  Vj^(t)r; 

6.  The  system  is  homogeneous  and  premixed; 

7.  The  energy  is  input  in  a  Gaussian  form  with  a  characeristic  radius 

R  which  increases  in  time;  and 
c 

8.  The  gas  never  gets  hot  enough  to  radiate  away  a  significant  fraction 
of  its  energy  during  the  Ignition  period. 


MODEL  COMPARISONS 


The  calibration  of  the  similarity  model  is  complete  when  the  Iccation 

at  which  the  Induction  parameter  should  be  evaluated  is  determined.  For 

this  purpose,  the  results  from  detailed  simulations  are  compared  to  those 

from  the  similarity  model.  Energy  deposition  in  both  models  is  linear  in 

time  at  a  rate  determined  by  requiring  an  energy  E^  to  be  deposited  in 

a  time  For  all  of  the  results  discussed  in  this  paper,  is  taken 

-4 

to  be  10  seconds.  As  in  the  similarity  model,  energy  deposition  in  the 
detailed  model  is  in  the  form  of  a  Gaussian  profile  in  space.  The  charac¬ 
teristic  radius  of  the  Gaussian,  R  ,  increases  with  time  and  is  determined  in 

c 

the  detailed  simulation  by  the  formula: 

2  ^  2/r^Ai(T/Ttt)dr 
^  J  ^  (T/Ta.)dr 


The  chemical  model  for  these  first  tests  was  taken  to  be  a  mixture 

of  H-:0-:N.  in  the  ratio  2:1:10  at  1  atm  and  T  »  300®K.  The  Induction 
Z  2  2  “ 

time  as  a  function  of  temperature  for  this  mixture  was  derived  from  the 
detailed  studies  of  the  ^2~^2  mechanisms  [13]  and  is  shown  in 

Fig.  2. 

In  the  similarity  model,  ic  was  estimated  by  comparing  the-'formula 
at  300°K, 

,S  8.4  X  10^ 

^M - 

o 

which  assumes  that  an  average  molecular  distance  a  and  an  average  mole¬ 
cular  weight  M  may  be  found,  to  the  more  exact  formulation 


M 


erg 


cm  sec  K, 


(14) 


»..i  i*i 


ITT  kjtj  “k  '^jk 


-1 


(15) 


where  is  a  function  of  {X^}  and  the  atomic  masses  {m^}  suggested  by 
Mason  and  Saxena  [ 19 ] .  This  gives  us 
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o  ••  3.16  A 


(16) 


M  -  2A.3. 

This  approximation  is  valid  because  the  similarity  solution  is  only  accurate 
for  ignition,  that  is,  before  any  major  amount  of  product  or  intermediates 
are  formed.  Then  the  parameter  k  is  found  using  the  definition  in  Eq.  (7). 

As  discussed  in  the  previous  section,  for  the  similarity  solution  k  is 
evaluated  at  the  central  (r»  0)  temperature.  In  contrast  the  detailed 
simulation  uses  Eq.  (IS)  to  evaluate  the  thermal  conductivity  at  each 
cell  at  each  time  step. 

In  the  first  case  studied,  the  initial  radius  of  deposition, 
uas  0.1  cm  and  both  models  were  configured  for  spherical  symmetry.  By 
varying  and  evaluating  the  induction  parameter  I(T,P^)  at  the  central 
temperature  (r=o)  the  similarity  solution  indicates  that  the  minimum  ignition 
energy  is  about  3.7mJ.  Figure  3  shows  the  t3rplcal  behavior  of  the  amplitude, 
A(t),  and  the  characteristic  radius,  R  ,  for  this  test  case.  If  I(T,P  )  is 
evaluated  at  R  ■  0.06  cm  the  minimum  ignition  energy  predicted  is  about  S.lmJ, 
and  if  it  is  evaluated  at  R^  ignition  is  not  predicted  even  when  E^  is  raised 
to  8mJ.  By  comparing  these  predictions  with  those  from  the  detailed  simula¬ 
tion  model  we  determine  the  location  at  which  I(T,P^)  must  be 
evaluated.  Results  from  the  detailed  simulation  are  shown  in  Fig.  4  for 
three  values  of  E^.  Ignition  occurs  when  is  greater  than  or  equal  to 
3.7mJ.  Therefore  good  agreement  between  the  two  models  for  the  case  under 
study  is  achieved  by  evaluating  the  induction  parameter  at  the  central 
temperature. 

The  central  temperature  and  the  Induction  parameter  have  been  shown 
for  three  values  of  E^  in  Figs.  5  and  6  respectively.  The  time  at  which 
ignition  occurs  (I  ■'1.0)  is  shown  by  a  in  Fig.  6.  A  comparison  with  the 
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results  from  the  detailed  simulation  (Fig.  4)  tells  us  that  both  models 
predict  Ignition  at  essentially  the  same  time  for  a  range  of  Input  energies. 

The  temperature  profiles  predicted  by  the  two  models  are  presented  In  Fig.  7. 

The  agreement  Is  very  good  except  for  long  times  when  the  detailed  simulation 
predicts  higher  temperatures  than  those  from  the  similarity  model.  This  Is 
primarily  due  to  the  heat  released  In  chemical  reactions  which  Is  not  Included 
In  the  similarity  model.  It  Is  Interesting  to  note  that  the  effects  of 
molecular  diffusion  included  In  the  detailed  simulation  model  are  not  con¬ 
sidered  In  the  similarity  solution  model.  The  results  presented  have  shown 
that  these  effects  are  not  as  Important  as  the  effects  of  thermal  conduction 
In  determining  the  Ignition  characteristics  of  the  system  under  Investigation. 

The  models  were  then  re-conflgured  for  cylindrical  geometry  and  the 
case,  R^=0,lcm,  was  again  investigated.  By  varying  E^,  the  similarity  solu¬ 
tion  Indicates  that  the  minimum  Ignition  energy  Is  about  3.5mJ/cm,  not  vejry 
different  from  the  spherical  case.  The  detailed  simulation  predicts  that 
the  minimum  Ignition  energy  is  between  3.3  and  3.7mJ/cm,  so  again  the  predictions 
of  both  models  are  comparable  when  the  Induction  parameter  is  evaluated  at 
r=0. 
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EFFECT  OF  QUENCHING  DISTAMCE 

As  discussed  in  the  previous  section,  good  agreement  between  the 
predictions  of  the  detailed  simulation  and  the  similarity  solution  is 
obtained  for  the  case,  »  0.1  cm.  In  order  to  determine  if  the  agree¬ 
ment  is  good  for  other  values  of  R  ,  various  cases  were  studied  and  these 

o 

have  been  summarized  in  Table  ll.  Evaluating  the  induction  peurameter  at  the 
central  temperature  provides  agreement  between  the  ignition  predictions 
of  the  two  models  for  all  the  cases  in  which  R^  is  greater  than  0.1  cm. 

However,  for  radii  smaller  than  0.1  cm  the  similarity  solution  predicts 
a  minimum  ignition  energy  that  is  lower  than  for  the  R^  =  0.1  case  described 
cd)ove.  In  fact,  as  the  initial  radius  of  deposition  decreases  below  0.1  cm 
the  similarity  solution  (with  l(T,Pg^  evaluated  at  a  fixed  radius)  predicts 
that  the  minimum  ignition  energy  also  decreases.  The  detailed  simulation 
does  not  agree  with  this  prediction.  The  radius  0.1  cm  appears  to  be  a 
"critical  radius"  below  which  the  minimvim  ignition  energy  again  increases. 

The  above  observation  is  similar  to  that  made  by  Blanc,  Guest,  von  Elbe  and 
Lewis  [20].  In  their  study  of  spark  ignition,  it  was  observed  that  there  was 
a  "critical  electrode  spacing"  below  which  the  minimum  ignition  energy  increased. 
The  critical  electrode  spacing  was  termed  the  "quenching  distour.r  a". 

To  study  the  phenomena  of  "absolute"  minimum  ignition  energies  and 
quenching  distances,  the  energy  deposition  in  the  detailed  simulation  needs 
to  be  modified.  Up  to  this  point  energy  was  deposited  in  a  Gaussian  profile 
with  a  characteristic  radius  which  increased  with  time.  This  was  necessary 
in  order  to  ensure  that  energy  was  deposited  in  a  manner  which  closely 
matches  the  one  derived  from  the  similarity  model.  In  the  similarity  model. 
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a  Gaussicui  profile  with  characteristic  radius  increasing  in  time  was 

essential  to  obtain  the  similarity  solution  described  earlier.  However, 

for  studying  the  effect  of  the  radius  of  deposition  on  the  minimum  ignition 

energy,  it  is  confusing  to  have  a  vaurying  radius  of  deposition.  Therefore, 

for  the  results  discussed  below,  energy  is  deposited  in  a  Gaussian  profile 

with  a  constant  radius  of  deposition.  As  before,  energy  deposition  is 

linear  in  time  at  a  rate  determined  by  requiring  an  energy  to  be  deposited 

in  a  time  t  . 

o 

For  a  particular  radius  of  deposition,  a  certain  amount  of  energy  is 
deposited  cuid  the  computations  are  carried  out  for  sufficient  time  until 
the  existence  (or  eibsence)  of  a  propagating  flame  is  definite.  By  repeating 
the  computations  for  different  values  of  E^  a  bound  for  the  minimum  ignition 
energy  for  that  particular  radius  is  obtained.  Similar  calculations  are 
performed  for  different  values  of  the  radius  of  deposition,  R^.  The  results 
of  such  investigations  are  shown  in  Pig.  8.  A  propagating  flame  results 
when  3.8  mJ  of  energy  is  deposited  in  a  sphere  with  a  radius  of  0.1  cm. 
However  if  the  same  amount  of  energy  is  deposited  in  a  sphere  of  smaller 
radius,  the  rate  of  heat  liberation  is  insufficient  to  compensate  for  the  rate 
of  heat  loss  and  consequently  there  is  no  ignition.  This  radius,  0.1  cm,  is 
the  "quench-radius"  for  this  particular  mixture.  For  radii  slightly  larger 
than  the  quench-radius,  the  minimum  ignition  energy  is  almost  constant  and 
for  larger  radii  (larger  than  0.11  cm)  the  minimum  ignition  energy  increases 
rapidly  with  increasing  radii.  Therefore  for  the  system  under  study,  the 
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absolute  minimviin  ignition  energy  is  about  3.7  mJ.  These  observations  are 

in  qualitative  agreement  with  those  of  Lewis  and  von  Elbe  [21].  Quantitative 

* 

comparisons  are  not  possible  since  the  composition  of  the  mixture  and  the  time 
for  energy  deposition  are  different. 

This  study  of  the  "quench-radius"  provides  the  insight  needed  to  explain 

the  discrepancy  in  the  predictions  of  the  similarity  solution  and  the  detailed 

model  for  small  radii.  We  conclude  that  a  vol\xne  smaller  than  the  "quench-volume" 

needs  to  be  maintained  at  a  temperature  T  for  a  time  which  is  longer  than  the 

corresponding  induction  time  for  ignition  to  occxir.  Therefore  the  concept 

of  an  induction  parameter  as  it  is  used  in  the  similarity  solution  is  not 

valid  for  very  small  radii.  If  the  absolute  minimum  ignition  energy  is  deposited 

in  a  sphere  of  radius  smaller  than  the  quench  radius,  ignition  will  not  occur 

since  the  heat  generation  rate  within  this  smaller  volume  cannot  compensate  for 

the  heat  loss  rate.  This  effect  is  not  significant  in  measurements  of  induction 

time  which  are  for  larger  radii.  Therefore  an  agreement  between  the  predictions 

of  the  two  models  can  be  forced  for  radii  smaller  than  the  quench-radius  by 

evaluating  the  induction  parameter  at  increasingly  larger  radii  (i.e.  lower 

temperatures).  For  example,  when  the  initial  radius  of  deposition  is  0.08  cm,  the 

detailed  simulation  predicts  that  the  minimum  ignition  energy  is  between  3.4 

and  5  mJ.  The  minimum  ignition  energy  predicted  by  the  similarity  solution 

is  3.4  mJ  if  the  induction  parameter  is  evaluated  at  0.06  cm  and  it  is 

5.1  mJ  if  the  induction  parameter  is  evaluated  at  0.08  cm.  For  the  case  when 

the  quench  radius,  r  »  0.1  cm,  the  formula 

^  BR 

r  »  r  -  Ae  (17) 

q 

A  =  2.84  X  10~^  and  B  =  58.78  gives  a  good  estimate  for  the  radius  at  which 
to  evaluate  the  induction  parameter.  We  are  currently  seeking  a  more  general 
formula  to  give  this  radius. 
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SUMMARY  AND  CONCLUSIONS 

This  paper  describes  a  theoretical  model  for  flar.;  ignition  based  on 

a  similarity  solution  which  may  be  used  to  predict  the  Ignition  properties 

of  a  homogeneous  mixture  of  gases.  The  model  requires  specification  of  the 

amount  of  energy  input,  E^,  the  length  of  time  over  which  the  energy  is 

deposited,  x  ,  and  the  radius  of  energy  deposition,  R  .  The  model  also 
o  o 

requires  basic  information  about  the  gas  mixture  which  includes  estimates 
of  the  thermal  conductivity  and  the  chemical  induction  time  of  the  gas 
mixture  as  functions  of  temperature.  Since  a  number  of  approximations 
have  been  applied  to  derive  the  model,  it  has  been  calibrated  and  its  range 
of  validity  determined.  This  has  been  done  by  comparing  its  predictions 
to  those  of  detailed  numerical  simulations. 

The  one-dimensional  detailed  numerical  simulation  used  solves  the 
set  of  coupled  partial  differential  equations  representing  conservation  of 
mass,  momentum,  and  energy  as  well  as  individual  species  densities.  For 
this  study  it  was  configured  with  an  open  boundary  at  one  end  to  simulate 
an  infinitely  large  system.  Energy  was  deposited  linearly  with  time,  as  it 
was  done  in  the  similarity  solution,  and  the  radius  of  deposition  was  chosen 
to  mimic  that  determined  by  the  characteristic  radius  in  the  similarity 
solution.  The  detailed  model,  however,  contained  calculations  of  the  thermal 
conductivity  and  chemical  kinetics  which  were  much  more  accurate  than  the 
approximations  used  in  the  similarity  solution.  Furthermore,  the  detailed 
model  contained  the  effects  of  molecular  diffusion,  which  were  not  at  all 
Included  in  the  similarity  solution. 

Predictions  from  the  similarity  model  consist  of  an  answer  as  to  whether 

or  not  a  given  mixture  would  ignice  given  E  ,  x  ,  and  R  .  Then  if  ignition 

o  o  o 

is  predicted,  the  model  gives  the  time  it  takes  for  ignition  to  occur.  In 
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contrast,  the  detailed  simulation  model  not  only  predicts  ignition, 
but  also  provides  the  structure  of  a  propagating  flame. 

Comparisons  of  the  predictions  of  the  two  models  show  that  the 
system  center  is  the  optimal  location  for  evaluation  of  an  induction 
parameter,  which  is  the  indicator  of  ignition  in  the  similarity  model. 
Fvurthermore ,  the  predictions  of  the  similarity  model  are  excellent  when 
energy  is  deposited  in  a  large  enough  voliame.  When  the  radius  of  deposition 
is  less  them  the  quench  radius  of  the  system,  the  similarity  model  is  not 
accurate:  the  similarity  model  predicts  ignition  when  the  detailed  model 
shows  it  does  not  occur.  This  can  be  compensated  for  by  evaluating  the 
induction  parameter  away  from  the  system  center.  However,  as  with  the 
value  of  the  quench  radius,  the  radius  at  which  the  induction  parameter 
must  be  evaluated  can  be  determined  by  comparisons  with  experiments  or 
calculations  performed  with  a  detailed  model. 

After  performing  a  series  of  studies  which  tested  and  calibrated  the 
similarity  model  in  cylindrical  and  spherical  geometries,  the  detailed 
model  was  reconfigured  to  evaluate  the  quench  radius  and  the  minimum  ignition 
energy  of  a  particular  gas  mixture.  This  study  provided  insight  into  the 
source  of  the  disagreement  between  the  predictions  of  the  similarity  and 
detailed  models  at  very  small  radii  of  energy  deposition.  The  error  occurred 
because  the  definition  of  an  induction  time  is  only  valid  for  a  volume  of 
material  that  is  large  and  homogeneous  enough  in  temperature  and  pressure 
so  that  diffusive  effects  are  not  important.  In  the  cases  where  the  radii 
of  deposition  were  very  small,  increased  thermal  conductivity  due  to  a  steep 
temperature  gradient  eroded  the  high  temperature  region  too  quickly.  Thus  in 
the  competition  between  chemical  energy  release  to  heat  a  region  and  thermal 
conduction  to  cool  it,  thermal  conduction  dominated. 


B-21 


Further  limits  of  validity  of  the  simplified  model  should  now  be  tested. 
These  include  investigating  the  sensitivity  of  predictions  to  time  or  mode 
of  energy  deposition,  to  variations  in  the  chemical  induction  time,  and  to 
the  effects  of  radiative  losses  discussed  in  the  Appendix.  These  kinds  of 
calibrations  and  determinations  of  the  limitations  zmd  sensitivity  of  the 
similarity  model  are  extremely  important  if  we  are  to  est^U^lish  the  level 
of  accuracy  required  for  input  data.  Once  this  has  been  done,  the  similarity 
model  can  be  used  with  experimental  input  data  to  estimate  quickly  whether 
a  material  is  flammable  and  if  so,  how  much  energy  is  required  to  ignite  it. 
Another  extension  of  the  work  presented  here  involves  using  the  detailed 
simulation  to  investigate  the  effect  of  geometry  and  method  of  energy 
deposition  on  ignition  energies  and  quenching  distances.  Both  these  studies 
are  currently  being  pursued. 
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APPENDIX 


The  (derivation  of  the  similarity  solution  model  is  discussed  in  detail 
here  for  the  case  of  spherically  symmetric  geometry.  Assuming  that  the 
pressure  field  is  essentially  constant  and  that  no  shocks  are  present,  the 
conservation  equations  for  momentum  emd  energy  may  be  combined  to  give 


dP  “It  f  IT 

^  »  0  =  -yPV  •  V  +  V  •  XNk„  7T  +  S(t)e  '  '  (Al) 

at  —  —  —  B  ““ 

^2  2 

whre  K  and  the  Gaussian  energy  deposition  term,  S(t)e  ,  have  already 
been  discussed.  The  other  equation  required  is  the  continuity  equation. 


1  ^ 
p  dt 


-7  •  V 


(A2) 


Assuming  dP/dt  *0,  Eq.  (Al)  gives  an  algebraic  equation  for  7^*  v  which 
may  be  combined  with  Eq.  (A2)  to  give 


T  dt  P  Y  _  *1.  _ 

CD 


{A3) 


where  P  is  the  background  pressure.  The  solution  is  then 

DD 


T(r,t) 

and 

P(r,t) 


(A4) 

(A5) 


where  T  and  p  are  the  background  temperature  and  density,  respectively. 

CD  CD 

Thus  the  nonlinear  slow  flow  equations  including  expansions  and  contractions 
of  the  flow  have  been  converted  into  a  single  equation  which  is  linear  in 
the  logarithm  of  the  temperature. 

The  total  energy  of  the  system  at  any  Instant  is  the  sum  of  the 
internal  energy  and  the  work  performed  in  expeuiding  the  heated  region.  It 
may  be  written 


1 

(A6) 


S  — I  ■  ■  ■  ■  ■ 

(V-l)k^(t) 

Y  ^flo  ' 

“  - r - P(A(t)).  (A7) 

(Y-1)  k'’(t) 

which  defines  the  integral  P{A(t)).  Differentiating  this  we  find  that 


E(t) 


P«Y 


Y-1 


f  [. . 


dE 

dt 


S(t) 

(Y-l)k^(t) 


(AS) 


which  may  be  equated  to 


dE  _  ^  at  5E  QA 
dt  ck  dt  dA  dt 


(A9) 


Thus  a  consistency  condition  has  been  specified  on  the  rates  of  change 
of  the  amplitude,  A(t),  and  the  scale  size  k"^^  (t)  for  the  heated  region. 

If  the  fluid  velocity  v  is  then  expcunded  such  that 

v(r)  »=^v^(t)  r,  (AlO) 

that  is,  only  the  linear  term  is  kept,  two  coupled  ordinary  differential 
equations  for  k  and  A  may  be  obtained. 


di-  “  ■  k^i  -2<k 


dt  Y  E 


(All) 


(A12) 


Then  the  expression  for  may  be  written  in  terms  of  the  integral  F(A) 


S 


P»(0)  -F*  (A)  ^ 


2  AF»(A)  -F(A) 


which  results  from  using  the  principle  of  energy  conservation  and  equating 
Eqs.  (A8)  and  (A9) . 

EFFECT  OF  RADIATIVE  LOSSES 

When  a  large  amount  of  energy  is  deposited  in  a  small  volume,  the  effect 
of  radiative  losses  may  be  significant.  Exact  calculation  of  all  the  low  tem¬ 
perature  radiation  effects  is  still  beyond  even  the  best  detailed  calculations 
today.  An  Idealized  model  can  be  Included  in  the  similarity  solution,  however, 
to  estimate  these  effects. 

Black  body  radiation  from  a  spherical  surface  at  radius  t  takes  energy 

out  of  the  system  at  a  rate 

dE  j  o  A 

— »  4irr  oT^(r).  (A14) 

dt 


In  the  similarity  solution  the  temperature  at  any  radius  r  is  given  by 

T  (r,  t)  -  T«  e*^*'^*  •  (A15) 


Therefore, 


dE 


rad 


dt 


,^42  4A(t)e 
4iroToo  re 


-k^(t)r^ 


(A16) 


However  at  t  “  0  (A  ••  0), 

dE  , 
rad 


dt 


-  0 


(A17) 


Therefore  the  background  radiation  into  the  volume  must  be  substracted  from 
(A  16)  to  give 


^  -  4.,Tt  ' 


/  2 
-  AiroTs  r 


(A18) 

The  next  step  is  to  choose  r  such  that  £q.  A18  is  maximized.  By  assuming 
that  the  loss  from  the  system  proceeds  at  the  maximum  rate  possible,  the  model 
is  sure  to  signal  the  onset  of  radiative  loss  effects  at  least  as  soon  as  they 
occur  in  the  physical  system  itself. 
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Table  I.  ®2”®2  Reactive  Mechanism 


-  AT®  exp  (-C/T)^*^ 


Reaction 

a(®> 

B 

c0>) 

Reference 

H  +  HO  0  +  H- 

1.40(-14) 

1.00 

3.50(+03) 

[241 

•*■  2 

3.00(-14) 

1.00 

4.48(403) 

[241 

H  +  H0_  i  H-  +  0- 

4.20(-U) 

0.00 

3.50(402) 

[241 

9.10(-11) 

0.00 

2.91(404) 

[241 

H  +  HO-  -*■  HO  +  HO 

4. 20 (-10) 

0.00 

9.50(402) 

[241 

2 

2.00(-U) 

0.00 

2.02(404) 

[24 1 

H  +  HO-  i  0  +  H-0 

8.30(-ll) 

0.00 

5.00(402) 

[251 

1.75(-12) 

0.45 

2.84(404) 

\ 

H  +  H-0-  i  HO  +  H 

2. 80 (-12) 

0.00 

1.90(403) 

[24l 

1.2a(-12) 

0.00 

9.40(403) 

[24] 

H  +  H-0-  i  HO  +  H-0 

5. 28 (-10) 

0.00 

4.50(403) 

[241 

2  2  2 

3.99(-10) 

0.00 

4.05(404) 

kt  - 

HO  +  H-  H  +  H-0 

1.83(-15) 

1.30 

1.84(403) 

[26l 

2  *■  2 

1.79 (-14) 

1.20 

9.61(403) 

[26l 

HO  +  HO  i  H-  +  0- 

1.09(-13) 

0.26 

1.47(404) 

k  -  kr/K 

•*■2  2 

2.82(-ll) 

0.00 

2.42(404) 

^  [27l  "" 

HO  +  HO  0  +  H-0 

1.00(-16) 

1.30 

0.00(400) 

[26] 

*■  2 

3.20(-15) 

1.16 

8.77(403) 

k^  -  k^/K^ 

HO  +  HO-  i  H-0  +  0- 

a. 30 (-11) 

0.00 

5.03(402) 

[28] 

2  -t-  2  2 

2.38(-10) 

0.17 

3.69(404) 

^ 

HO  +  H-0  i  HO  +  H 

1.70 (-11) 

0.00 

9.10(402) 

[24] 

2  2  2 

4. 700-11) 

0.00 

1.65(404) 

[24l 

HO  +  0-  i  HO-  +  0- 

1.60 (-12) 

0.00 

9.56(402) 

[25l 

3  2  2 

6. 69 (-14) 

0.33 

2.04(404) 

\  /K^ 

HO  +  H-  i  HO  +  H-0 

1.20 (-12) 

0.00 

9.41(403) 

[27] 

2  2 

1.33(-14) 

0.43 

3.62(404) 

\  -  V*^C 

HO,  +  HO-  -*•  H-0-  +  0- 

3.00(-ll) 

0.00 

5.00(402) 

[25] 

2  2^22  2 

1.57 (-09) 

-0.38 

2.20(404) 

k  -  k-/K 
r  £  c 
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Table  I.  (continued) 

H--0.  Elementary  Reactive  Mechanism 

A  4m 


Reaction 

AT  exp  1 

B 

C-C/T)^®^  X 

Reference 

0  +  HO  -*■  H  +  0- 

2. 72 (-12) 

0.28 

-8.10(+01) 

K  - 

3.70(-10) 

0.00 

8.45 (+03) 

X  .  r  c 
(241 

0  +  HO  J  HO  +  0- 

8. 32 (-11) 

0.00 

5,03(+02) 

[28] 

2, 20 (-11) 

0.18 

2. 82  (+04) 

k  -  k,/K 
r  f  c 

0  +  H,0,  t  H^O  +  0« 

1.40 (-12) 

0.00 

2. 12 (+03) 

[25] 

t 

5. 70 (-14) 

0.52 

4.48(+04) 

kr  -  k,/K^ 

0  +  H_0,  2:  HO  +  HO, 

1.40 (-12) 

0.00 

2.13(-M)3) 

[25 1 

2. 07 (-15) 

0.64 

8.23(-H)3) 

-  k-/K 
r  X  c 

H  +  H  +  MiH,  +  M 

1.80 (-30) 

-1.00 

O.OO(-H)O) 

[241 

3.70(-10) 

0.00 

4. 83 (-04) 

[241 

H  +  HO  +  M  i  H-0  +  M 

6. 20 (-26) 

-2.00 

0.00 (+00) 

[24] 

5. 80 (-09) 

0.00 

5.29(+04) 

[24] 

H  +  0,  +  M  2:  HO,  +  M 

4. 14 (-33) 

0.00 

-5.00(-H)2) 

[24  J 

3. 50 (-09) 

0.00 

2.30(+O4) 

P4] 

HO  +  HO  +  M  i  H,0,  +  M 

2. 50 (-33) 

0.00 

-2.55(-M)3) 

[24] 

2. 00 (-07) 

0.00 

2. 29 (+04) 

[24] 

o  +  h  +  m2:ho  +  m 

8. 28 (-29) 

-1.00 

0.00(-f00) 

[29] 

2. 33 (-10) 

0.21 

5. 10  (104) 

k  -  k,/K 
r  f  c 

0  +  HO  +  M  HO,  +  M 

2. 80 (-31) 

0.00 

0.00(100) 

[29) 

1.10 (-04) 

-0.43 

3.22(104) 

k  -  k-/K 
r  £  c 

O  +  O+  MiO,  +  M 

5.20(-35) 

0.00 

-9.00(102) 

[24] 

3. 00 (-06) 

-1.00 

5.94(104) 

[24] 

(a)  Bimolecular  reaction  rate  constants  are  given  In  units  of  cm  /(molecule  sec). 


6  2 

Tezmolecular  reaction  rate  constants  are  given  In  units  of  cm  /(molecule  sec), 
(b)  Exponentials  to  the  base  10  are  given  in  parenthesis;  l.e.,  I.OO(-IO)  - 
1.00  X  10”^°. 
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Table  IZ 


Case 

ZnitlaJ.  Radius 

of  Deposition 

R  (cm) 
o 

Total  Energy 

Deposited 

E  (na) 
o 

Ignition  Prediction 
Detailed  Similarity 

Numerical  Solution 
•  Simulation 

1 

0.1 

4.0 

yes 

yes 

3.7 

yes 

yes 

3.0 

no 

no 

2 

0.11 

5.0 

yes 

yes 

3.0 

no 

no 

3 

0.12 

7.0 

yes 

yes 

6.0 

no 

no 

4 

0.09 

3.3 

no 

yes 

1 

5 

0.08 

3.3 

no 

j 

yes  1 
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FIGURE  CAPTIONS 


Figure  1 


Figure  2 


Figure  3 


Figure  4 


Figure  5 


Figure  6 


Figure  7 


Figure  8 


Time  history  of  the  temperature  profile  in  a 

2:1:10  mixture  as  predicted  by  the  detailed  numerical  simulation 
model. 

Temperature  as  a  function  of  induction  time  for  H2:02:N2/2:1:10 

mixture  evaluated  using  the  chemical  reaction  scheme  in  Table  I. 

The  nonlinear  amplitude.  A,  2md  the  characteristic  radius,  R  , 

c 

as  functions  of  time  calculated  using  the  similarity  solution 
model. 

Calculations  of  the  central  teii5>erature  as  a  function  of  time 
for  three  values  of  using  the  detailed  numerical  simluation 
model. 

The  central  temperature  as  a  function  of  time  for  three  values 
of  E^  as  predicted  by  the  similarity  solution  model. 

Calculations  of  the  induction  parameter  as  a  function  of  time 
for  three  values  of  E^  using  the  similarity  solution  model. 

The  indicates  the  predicted  time  of  ignition. 

Comparisons  of  the  time  history  of  the  temperature  profiles 
predicted  by  the  detailed  numerical  model  and  the  similarity 
solution  model. 

The  minimum  ignition  energy  as  a  function  of  the  radius  of 
energy  deposition  calculated  using  detailed  numerical  simulations. 
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